package com.sr.srinterpolation.util;

import android.util.Log;

/**
 * 拉格朗日插值法
 */
public class LagrangeInterpolationUtil {
    private static final String TAG = "LagrangeInterpolationUtil";

    /**
     * 返回拉格朗日插值多项式
     *
     * @param x {1, 2, 3, 4, 5, ... , n}
     * @param y {1, 2, 4, 8, 16, ... , 114514}
     * @return {{分子数组}{分母数组}}
     */
    public static int[][] polynomial(int[] x, int[] y) {
        // 计算拉格朗日插值多项式每项的分母部分
        int[] lDens = lDens(x);
        // 计算拉格朗日插值多项式通分后的分母
        int lReduceDen = Math.abs(lcm(lDens));
        // 计算拉格朗日插值多项式每项的通分倍率
        int[] lRadios = reduceRadio(lDens, lReduceDen);
        // 计算分子部分
        int[] nums = nums(x, y, lRadios);
        // 约分得出最终结果
        int[][] result = MathUtil.reductionOfAFraction(nums, lReduceDen);

        for (int num : nums) {
            Log.i(TAG, "num: " + num);
        }
        for (int resultNum : result[0]) {
            Log.i(TAG, "resultNum: " + resultNum);
        }
        for (int resultDen : result[1]) {
            Log.i(TAG, "resultDen: " + resultDen);
        }
        return result;
    }

    // 计算分子部分
    // tips: n个数的拉格朗日插值多项式的第一项x部分如下：
    // 3: (x - x1)(x - x2) = x^2 - (x1 + x2)x + x1x2
    // 4: (x - x1)(x - x2)(x - x3) = x^3 - (x1 + x2 + x3)x^2 + (x1x2 + x1x3 + x2x3)x - x1x2x3
    // 5: (x - x1)(x - x2)(x - x3)(x - x4) = x^4 - (x1 + x2 + x3 + x4)x^3 + (x1x2 + x1x3 + x1x4 + x2x3 + x2x4 + x3x4)x^2 - (x1x2x3 + x1x2x4 + x1x3x4 + x2x3x4)x + x1x2x3x4
    private static int[] nums(int[] x, int[] y, int[] lRadios) {
        int n = x.length;
        int[] nums = new int[n];
        // 以合并同类后的多项式中每项的系数为循环，i为该项的次数
        for (int i = 0; i < n; i++) {
            int num = 0;
            // 该次项中每项中乘数的个数
            int count = n - i - 1;
            // 以i次项中每个源自拉格朗日插值多项式中的分子为循环，j为拉格朗日插值多项式中的第几项
            for (int j = 0; j < n; j++) {
                int lNum = subCoefficient(x, j, 0, 1, count);
                num += lNum * y[j] * lRadios[j];
            }
            // 将i次项系数分子存入数组
            nums[i] = num;
        }
        return nums;
    }

    /**
     * 计算分子中拉格朗日插值多项式每项中i次项的系数
     *
     * @param x         x数组
     * @param lIndex    拉格朗日插值多项式的项数
     * @param index     遍历数组的位置
     * @param sub       系数的每个乘数项（初始传入1）
     * @param needDepth 剩余层数
     * @return 分子中拉格朗日插值多项式每项中i次项的系数
     */
    private static int subCoefficient(int[] x, int lIndex, int index, int sub, int needDepth) {
        if (needDepth == 0) {
            return sub;
        }
        int result = 0;
        for (int i = index; i < x.length - needDepth + 1; i++) {
            if (i == lIndex) continue;
            result += subCoefficient(x, lIndex, i + 1, sub * -x[i], needDepth - 1);
        }
        return result;
    }

    // 计算拉格朗日插值多项式每项的分母部分
    private static int[] lDens(int[] x) {
        int n = x.length;
        int[] dens = new int[n];
        for (int i = 0; i < n; i++) {
            int den = 1;
            for (int j = 0; j < n; j++) {
                if (i != j) {
                    den *= x[i] - x[j];
                }
            }
            dens[i] = den;
        }
        return dens;
    }

    // 获取通分倍率（所有分数通分后，每项所需乘的倍数）
    private static int[] reduceRadio(int[] dens, int reduceDen) {
        int n = dens.length;
        int[] radios = new int[n];
        for (int i = 0; i < n; i++) {
            radios[i] = reduceDen / dens[i];
        }
        return radios;
    }

    // 多个数的最小公倍数
    private static int lcm(int[] a) {
        int lcm = MathUtil.lcm(a[0], a[1]);
        for (int i = 2; i < a.length; i++) {
            lcm = MathUtil.lcm(lcm, a[i]);
        }
        return lcm;
    }
}
